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The increasing penetration of electric vehicles over the coming decades, taken together 
with the high cost to upgrade local distribution networks and consumer demand for home 
charging, suggest that managing congestion on low voltage networks will be a crucial com¬ 
ponent of the electric vehicle revolution and the move away from fossil fuels in transporta¬ 
tion. Here, we model the max-flow and proportional fairness protocols for the control of 
congestion caused by a fleet of vehicles charging on two real-world distribution networks. 
We show that the system undergoes a continuous phase transition to a congested state as 
a function of the rate of vehicles plugging to the network to charge. We focus on the or¬ 
der parameter and its fluctuations close to the phase transition, and show that the critical 
point depends on the choice of congestion protocol. Finally, we analyse the inequality in 
the charging times as the vehicle arrival rate increases, and show that charging times are 
considerably more equitable in proportional fairness than in max-flow. 
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I. INTRODUCTION 


Electric vehicles may become competitive, in terms of total ownership costs, with internal- 
combustion engine vehicles over the next couple of decades. Studies in the United States and the 
UK suggest the current power grid has enough generation capacity to charge 70% of cars and light 
trucks overnight, during periods of low demand tit]. A recent survey suggests, however, vehicle 
owners prefer home charging, would consider charging their vehicles during the day (typically 
between 6 and 10 pm), and are unwilling to accept a charging time of 8 hours t2D. The time to 
fully charge the battery of an electric vehicle at home currently varies from 18 hours (Level 1, in the 
United States at 110 V and 15 A with a charge power of 1.4 kW) to 4 hours (Level 2, at 220 V, 30 
A with a charge power of 6.6 kW). Alternatively, electric vehicles could charge at public charging 
stations at Level 3 in less than 30 minutes yD. Taken together, consumer behaviour and advances 
in battery technology may lead to a rise in peak demand with the increasing penetration of electric 
vehicles, overloading distribution networks and requiring local infrastructure reinforcement lO- 
lll]. To reduce the cost of upgrades to the last mile of cables, network operators may need to 
coordinate charging strategies in a way that is both technically and socially acceptable. To achieve 
this goal, network designers could implement charging protocols that prioritise the access of a 
fleet of electric vehicles to electric power, thus simultaneously managing network congestion and 
accounting for the fairness of user allocations. 

Through a series of papers, the power grid has recently gained increased visibility in the sci¬ 
entific community [8, and physicists have helped to increase our understanding of its synchro¬ 


nization llld 


transitions 


0 


Land stability |ll2L 


13h . In parallel, recent advances in optimization and phase 


15h suggest that the tools of critical phenomena and optimization can be merged. 


opening up new horizons. From the point of view of the distribution network operator, the problem 
of vehicle charging is to manage congestion on distribution networks, while respecting Kirchholf’s 
laws and keeping voltage drops bounded. Here, we explore two congestion control mechanisms: 
max-flow and proportional fairness. We show that if too many vehicles plug-in to the network, 
charging takes too long, more cars arrive than leave fully charged, and the system undergoes a 


continuous phase transition to a congested state lla 


ITy, where the critical point depends on the 


choice of congestion control algorithm. By gaining insights into the critical behaviour that nat¬ 
urally emerges with the increase of the number of vehicles, we hope to help network designers 
decide which algorithms to implement in the real-world. 
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II. THE MODEL 


Physicists are familiar with simulated annealing, a global optimization method that can avoid 
becoming trapped in a local optimum. In principle, it converges to the global optimum, but in prac¬ 
tice this is not guaranteed (see e.g. Ill 8142111 1 because the required theoretical cooling schedules are 
too slow to use in implementations. In contrast, convex optimization always finds the solution, 
if it exists, independently of the starting point. Convex optimization problems can be solved effi¬ 
ciently (typically in polynomial time), even for problems with hundreds of variables and thousands 
of constraints, using interior-point methods 1I22I1 . The burgeoning field of convex optimization in 
electricity networks 112314250 is a good example of an application of the mathematical framework 
developed over the last 20 years. Indeed, the extensive numerical simulations we present here are 
only possible due to techniques developed since 2012 ll24l426ll . The networks that we study are 
relatively small. The stochasticity of vehicle arrival times, however, implies solving an optimiza¬ 
tion problem in each time step if the state of the system changes. Hence, to gain insights into the 
steady state of vehicle charging, efficient algorithms are a necessity at the design stage. Of course, 
real-world implementations also depend on efficient algorithms, which would need to run online, 
often in large urban distribution networks. 


An optimization problem is determined by a function of a set of variables (the objective func¬ 
tion), for which we seek a minimurn, and a set of upper bound constraints that restrict the domain 
{or feasible set) of those variables 1I22I1 . A point is feasible if it belongs to the feasible set, and is 
optimal if it is the minimum of the objective function in the feasible set. An optimization problem 
is convex if both the objective function and the constraints are convex, in which case the objective 
function has a global minimum. A convex relaxation of an optimization problem P is a convex 
optimization problem P' with an enlarged feasible set. If the optimum of P' is feasible for P, 
it is also the optimum for P and we say the relaxation is exact. Hence, convex relaxations are 
more attractive than approximate methods, such as linearisations, because the feasibility of the 
relaxed optimum of P', which can be verified either analytically or numerically, is a certificate of 
the exactness of the relaxation. 


Consider a tree topology, such that electric power is distributed from a root node to electric 
vehicles that charge at the nodes. Let P{t) be the feasible set of power allocations at time t, i.e. the 
set of all allocations of power to electric vehicles that do not violate the operational constraints 
of the distribution network. Each feasible allocation P(t) e P{t) is a vector P(t) = {Pi{t) : / = 
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1,.. where N(t) is the number of vehicles in the network at time t. Vehicle / derives 

a utility Ui(Pi(t)) from the allocated charging power Pi(t), and we wish to select the allocation 
that maximises the sum of vehicle utilities 1270. This allocation acts as a network protocol that 
distributes network capacity among users, and solves the following problem: 


N{t) 


maximise E Ui{P,(t)) 


(la) 


/=i 


subject to P{t)eP{t). (lb) 

Here we explore two user utility functions. First, we consider the non-unique max-flow allocations 
given by Ui iPi(t)) = Pi(t), i.e. we maximise the instantaneous aggregate power sent from the root 
node to the vehicles, which is a benchmark of efficient network throughput feSfl . Such allocations, 
however, can also leave users with zero power, which is considered unfair from the user point of 
view. Hence, we next consider the proportional fairness allocation. Mathematically, the problem 
is to find the feasible allocation that maximises the sum of the logarithm of user rates, that is 
Ui {Piif)) = logfPiif)). The proportional fairness allocation is especial, because the users and the 
network operator simultaneously maximise their utility functions s. Furthermore, the problem 
is convex, and so can be solved in polynomial time 1220, and it can be naturally extended by adding 
positive weights to each term in the objective function Eq. (fTah . to account for diversity in user 
demand or for more than one user at some nodes KTfl . For the compact and convex set Pit), it can 
be shown that the allocation P^^{t) that maximises Eq. (fTal) . satisfies 11271 

Pi{t)-P^^it) 


290: 


z 


< 0 . 


( 2 ) 


/=i 

This allocation is known as proportionally fair, because the aggregate of proportional changes 
with respect to all other feasible allocations is non-negative. In other words, Eq. Q implies that 
to increase the instantaneous power allocated to a vehicle by a percentage e, we have to decrease a 
set of other power allocations, such that the sum of the percentage decreases is larger or equal to e. 
In contrast, in max-flow, to increase the instantaneous power allocated to a vehicle by e, we have 
to decrease the sum of instantaneous powers allocated to other vehicles at least by e. It turns out 
that fairness and congestion control are two sides of the same coin, because the existing algorithms 
for fair allocations also manage network congestion 11271 boUs511. Moreover, in the analysis of the 
parallel problem for communication networks, proportional fairness has emerged as a compromise 
between efficiency and fairness with an attractive interpretation in terms of shadow prices and a 
market clearing equilibrium 11271 Section 7.2]. 
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FIG. 1. Schematic illustration of (a) a distribution network, (b) the circuit of a network edge. Electric 
vehicles choose a charging node with uniform probability, and plug-in to the node until fully charged, as 
illustrated by the electric vehicle icons on the network. Network edge eij has impedance Zjj - Rjj + iXij. 
The power consumed by the subtree iti (y) rooted at node j (area shaded in purple) is S fti(y) = F’fti(y) -i- iQ(t\{j), 
where vehicles consume real power only, but network edges have both active (real) and reactive (imaginary) 
power losses. 


The simplest flow model in electricity networks is the DC power flow, which is a linearisation 
of the AC power flow equations, and thus can be solved with tools of linear programming. It 
assumes that voltage amplitude is constant on all nodes, a good approximation at the level of the 
high-voltage transmission network, but a poor one on local distribution networks. Indeed, voltage 
drops are significant in distribution networks, and determine the network capacity, which leads us 
to using models of power flow specific to distribution networks 11361] . We abstract the distribution 
network to a rooted directed tree rti (r) with node (often called bu^ set 'V, edge (also called branch) 
set fi, and a root node r (feeder) that injects power into the tree 111. Edge e 6 connects node i to 
node j, where i is closer to the root than j, and is characterised by the impedance Z;y = Rij + iXjj, 
where is the edge resistance and Xij the edge reactance. The power loss along edge is 
given by Sift) = Pift) -i- iQift), where Pi ft) is the real power loss, and Qift) the reactive power 
loss. Electric vehicle / receives active power Pft) until charged, but does not consume reactive 
power —see Fig.dl^a). The vector V(t) denotes the voltage allocated to the nodes. The voltage 


drop AVij down the edge is the difference between the amplitude of the voltage phasors V, and 


* We write iti (r) instead of fti {'V, fi, r) to simplify the notation. 
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Vj, assuming node i is closer to the root r than node j |^. Let iti (j) denote the subtree of the 
distribution network rooted in node j, with node set "Vihij) and edge set firti(y). Let P^(j) denote the 
active power, and Q^(j) the reactive power consumed by the subtree rti (j )—see Fig.lH Kirchhoff’s 
voltage law applied to the circuit in Fig.[Hb) yields (see Appendix lAh: 


- V](t) - = 0. 


(3) 


Vehicle / has a battery with capacity B that charges with the instantaneous power Pi(t) from 
empty (at arrival time) to full (at departure time), and the level of battery charge is the time integral 
of instantaneous power. Vehicles arrive to the network, choose a node to charge randomly with 
uniform probability, charge until their battery is full, an lastly leave the network. At each time step, 
the network solves the congestion control problem to allocate instantaneous power to the vehicles. 
The max-flow problem maximises the instantaneous aggregate power sent from the root node to 
the electric vehicles, respecting the constraints of distribution networks: the voltage drop along 
edges obeys Eq. ©, and node voltages are within ((1 - a)V^mirmi, (1 + aWnominui) for a e (0,1), 
with a = 0.1 typically 11361] . Thus, to find the max-flow allocation of power to the vehicles, we 
solve the optimization problem for fixed t: 


N(t) 


maximise U{t) = ^^Pi(t) 

i=\ 


(4a) 


subject to (1 Cp)l^nominal — F,(t) ^ (1 Ct^Vnominah ^ ^ *3^ (4b) 

Vi(tWj(t) - Vj(tWj(t) - PHJ)mij - QHJ)mij = 0, e 6. (4c) 

Constraint (l4bh sets the limits on the node voltage. Equation (l43 is the physical law coupling 
voltage to power, generalized from Eq. © for the subtree rti (j), and encodes both Kirchhoff’s 
voltage law on network edges and Kirchhoff’s current law applied recursively at each node of the 
subtree (see Appendix |B]). We do not need to apply Kirchhoff’s voltage law on network loops, 
however, because local distribution networks are approximately trees, and thus are cycle-free. 
Constraint (l4^ is quadratic, hence not convex in general, which implies that the problem is not 
directly solvable by polynomial time methods. To overcome this limitation, we make a change 
of variables in problem © by defining a weighted adjacency matrix W(t), such that edge eij 


corresponds to the 2x2 principal submatrix W(eij, t) defined by 11381 


m 






Vf(t) 



'WuiO Wijit) 

VjitWiit) 

vj(t) ^ 


.Wjiit) Wjjitx 


( 5 ) 
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where Wij(t) = Wpit), because Vj(t) e R. The matrices W(eij, t) are positive semidefinite, 


because their eigenvalues (Ti = 0 and A 2 = + V]) are non-negative, and rank one because 


they are of the form vv^. Hence, constraint (1^ can be replaced by three constraints: the first 
substitutes the quadratic terms in the voltages with linear terms in the t), and the second and 
third constraints guarantee that the W{eij, t) are positive semidefinite and rank one. 

The solution of problem Q is on the Pareto frontiei^ since we maximise an increasing func¬ 
tion in the obiective. The rank one constraint is nonconvex, but it does not change the Pareto 

nn 

frontier or the optimum ll25L 140(1 . and we remove it to relax problem Q to: 


m) 


maximise 

Wit) 

u(t) = YjPi(t) 

i=i 


( 6 a) 

subject to 

((1 ~ tr)Vfiominal) — — ((1 "i” tr)Vnominal) » 


( 6 b) 


Wij(t) - - P<,u)(t)Rij - = 0 , 

Bij e 6 

( 6 c) 


W(eij, t) > 0 , 

^ij S fi. 

( 6 d) 


where the generalized inequality in constraint (IMI) means the W{eij,t) matrices are positive 
semidefinite 

The problem of allocating power to vehicles in a proportional fair way has the same constraints 
as problem ®, however, the objective function is the sum of the logarithm of the power. It turns 
out, however, that it is computationally more efficient to aggregate vehicles at the nodes, and to 
maximise the sum of power allocated to the nodes, rather than the vehicles. To show this, we 
observe that all vehicles are equivalent, and thus the power Pi (0 allocated to node i is divided 
equally among the vehicles charging on the node at each time step. Hence, if one or more vehicles 
is charging on node i, each gets the instantaneous power: 


Pi{t) = 


Pi(0 

W,(0’ 


(7) 


where w ,(0 = ^//(O is the number of electric vehicles charging on node i at time t, and 

A, 7(0 = 1 if electric vehicle / is charging on node i at time t and zero otherwise. Hence, the 


^ We say that a power allocation {Pi} for / = 1,..., is better than another if Pi > P'j for all I, and for some I, 
Pi > P'l- A power allocation is Pareto optimal or efficient if there is no better power allocation. The Pareto frontier 
of a set is the set of all Pareto optimal points. 
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proportional fair allocation is given by (see Appendix |Cl): 


maximise 

W(t) 

U(0= ^ w,(0iogPi(0 

ie^+ 



(8a) 

subject to 

((1 U() ^nominal) — ^ iiif) — ((1 "1" ^)^nominal) 5 


(8b) 


- Q^mij 

= 0, 

Cij e 6 

(8c) 


WiCij, t) > 0, 


^ij S fi. 

(8d) 


where is the subset of nodes with at least one vehicle, and we can recover the instantaneous 
power allocated to electric vehicle /, located at node i, from Eq. ©. The complexity of the prob¬ 
lem ((SI) thus scales with the number I'T'I of nodes, which is typically smaller than the number N(t) 
of vehicles for large arrival rates A. Similarly, we also aggregated vehicles in the implementation 
of problem db]), but omit the proof. 

To study the behaviour of max-flow and proportional fairness as a function of the number of 
vehicles arriving at the network to be charged, we implement a discrete simulator that solves 
the congestion control problem in discrete time steps, starting with no vehicles charging on the 
network. Vehicles arrive at the network in continuous time (following a Poisson process with rate 
T) and with empty batteries, choose a node with uniform probability amongst all nodes (excluding 
the root), and charge at that node until their battery is full, at which point in time they leave the 
network. Once a vehicle plugs into a node, the congestion control algorithm will allocate it an 
instantaneous power, which is a function of the network topology and electrical elements, as well 
as the location of other vehicles. 


At each time step, we first check whether the number of charging vehicles changed (i.e. vehicles 
left the network fully charged, or new vehicles arrived to be charged), and if it has, we solve the 
max-flow problem ((3) and the proportional fairness problem ([8]), which allocate a constant power 
during the time step to each of the charging vehicles. Next, we update the status of batteries at 
the end of the time step. The simulation terminates when the simulation time reaches the time 
horizon. We simulated vehicles charging on the realistic SCE 47-bus and SCE 56-bus distribution 
networks ll38l] . which are detailed in Fig. |2l To characterise the system behaviour in detail, we 
varied the arrival rate A from 0 to 1 in steps of 0.05 (0.005 close to the critical points), and for 
each A value we simulated an ensemble of 25 independent realisations of simulation runs, each 
simulation running for 15,000 time units (150,000 time units close to the critical point). We ran 
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FIG. 2. Topology of the (a) SCE 47-bus and (b) S CE 56-bus networks. Node indexes identify the edges, and 


edge resistance and reactance is taken from 1138h . Node 1 is the root node in both networks. Nodes 13, 17, 


19, 23 and 24 of the SCE 47-bus network (in lighter colour) are photovoltaic generators, and we removed 
them from the network. 


the simulations using CVXOPT Il42n on the ETHZ Brutus elusterlj due to the high eomputational 


0 ( 


requirements. The computational time is comparable for max-flow and proportional fairness and 
for the 47-bus and 56-bus networks, but it is grows with A. For example, to simulate 5,000 time 
units of the proportional fairness algorithm for A = 1.0 on the 47-bus network takes approximately 
40 hours, but 4 minutes for A = 0.05 . 

We set the battery capacity 5 = 1 for all vehicles, and the nominal voltage Vnominai = 1- Scaling 
ynominal by jA, for yS G (0, oo), implics scaling both the the power delivered to vehicles and the 
battery capacity by To see this, observe that problems © and dS]) are invariant upon the scaling 
"^'nominal = fiVnominai, V'(t) = /3Vi(t) for all nodal voltagcs, P'jit) = p^Piit), P'^it) = 0^P(nit) and 


Q'^{t) = P and B' = /A B. Considering these scaling properties, our simulations can be 

extended to values of Vnominai ^ C provided the vehicle capacity B is rescaled accordingly, and we 
use this property to rescale the problem when convenient. 


III. NUMERICAL RESULTS 


We find critical behaviour that resembles results found in communication networks, in that both 


systems undergo a continuous phase transition [1430. In order to characterize this phase transition 


we adopt the order parameter ri{A) that represents the ratio at the steady state between the number 


^ https://wwwl.ethz.ch/id/services/list/comp_zentral/cluster/index_EN 
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FIG. 3. (a) Order parameter 77 as a function of the arrival rate A, for the SCE 56-bus (filled symbols) and 47- 
bus (unfilled symbols) networks, where we apply the max-flow (circular symbols) and proportional fairness 
(square symbols) algorithms for the simulation horizon of 1.5x 10"^ time units. We plot a zoom of the critical 
region for the (b) 56-bus network and (c) 47-bus network for the longer horizon of 10^ time units. Panel 
(c) suggests the critical arrival rate is different for the max-flow and proportional fairness algorithms in the 
47-bus network. Symbols show average values over an ensemble of 25 runs and shaded areas represent 95% 
confidence intervals. 


of uncharged vehicles and the number of vehicles that arrive at the network to be charged [1430: 

1<AA^(0) 

; 7 (i) = lim-———, (9) 

A At 


where AN(t) = N(t + At) - N(t) and (...) indicates an average over time windows of width At. 
We calculate t](A) in the steady state, that is limt^co AN(t) oc At. For arrival rates A < A^, all 
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vehicles that plug-in to the network with empty batteries within a large enough time window leave 
fully charged within that period (free flow phase), but for A > Ac some vehicles have to wait for 
increasingly long times to fully charge (congested phase). The order parameter characterises the 
phase transition: ri(A) = 0 in the free-flow regime, and ri(A) > 0 in the congested phase, a higher 
order parameter meaning that queues of charging vehicles build up more rapidly. 

Figure [3] is a plot of the order parameter for the 47 and 56-bus networks and the two conges¬ 
tion control methods, as a function of A. Simulation results shown in Figs. |3]^a) suggest that Ac 
depends on several factors (the network topology, the complex impedance on the edges, battery 
capacity, Vnominah as well as the position of vehicles on the network). At this resolution of the 
control parameter, it is unclear, however, whether the critical point is the same for max-flow and 
proportional fairness in both networks. To clarify this, we studied the order parameter with higher 
resolution close to the critical points—see Figs. I3b) and (c). The critical point is numerically 
indistinguishable for max-flow and proportional fairness in the 56-bus network. In the 47-bus 
network, however, we find that Ac is larger for proportional fairness than for max-flow. 

The number N(t) of charging vehicles at time t fluctuates widely close to the critical point, 
and thus it is difficult to determine Ac from Fig. [3l To overcome this limitation, we adopt the 
susceptibility-like function MSh : 


X(A) = lim At cr^(At), (10) 

At^oo 

where At is the length of a time window, and crj^(At) is the standard deviation of the order parameter 
rj. To compute we first consider a long time series and split it into windows with length 
At. We next determine the value of the order parameter in each window, and finally calculate 
the standard deviation of these values. The susceptibility displays a singular point at Ac (see 
Fig.® , allowing us to study the dependencies of the critical arrival rates on the congestion control 
algorithm, as well as network topology and size. 

Similarly to our analysis of ri(A), the values of Ac are indistinguishable in the 56-bus network. 
In contrast, however, in the 47-bus network the singular point of;t'(/l) is smaller for max-flow than 
for proportional fairness. This suggests that proportional fairness charges a slightly larger number 
of vehicles than max-flow, and is thus marginally more efficient, on a neighbourhood of its critical 
point. To support this conclusion, we show in Figs. |4tc) and (d) four representative instances of 
the time series of the number of vehicles charging on the 47-bus network at T = 0.39. The number 
N(t) of vehicles grows linearly with time in max-flow in all four cases, suggesting that the critical 
point is below A = 0.39 for this algorithm. In contrast, N(t) oscillates in proportional fairness. 
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FIG. 4. Susceptibility as a function of the arrival rate A, for the for the (a) SCE 56-bus (filled symbols) 
and (b) 47-bus (unfilled symbols) networks, where we apply the max-flow (circular symbols) and propor¬ 
tional fairness (square symbols) algorithms for the time horizon of 10^ time units. Vertical lines show the 
value of the critical points for max-flow (ME) and proportional fairness (PE). Panel (b) shows the differ¬ 
ence in the critical arrival rate for the two congestion control algorithms. To illustrate this difference, we 
plot in (c) and (d) representative time series for the 47-bus network for A = 0.39, showing that, within the 
time horizon, max-flow is supercritical, whereas proportional fairness is subcritical. Symbols show average 
values over an ensemble of 25 runs and shaded areas represent 95% confidence intervals. 


suggesting that the critical point is above A = 0.39, in agreement with the analysis ofx('i)- 

The two congestion control algorithms lead to different allocations of instantaneous power, with 
vehicles charging in different order and over different time intervals. If there are vehicles on a path 
p between the root and a leaf node, the voltage drops with increasing distance from the root, the 
lower limit voltage constraint (1^ is fulfilled at equality for one node on p, and nodes further away 
than that will not receive any power. The objective function of proportional fairness guarantees 
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FIG. 5. Gini coefficient G of the charging time as a function of the electric vehicle arrival rate A, for the 
SCE47-bus (unfilled symbols) and 56-bus (filled symbols) networks, where we apply the max-flow (circular 
symbols) and proportional fairness (square symbols) algorithms. We run the simulation for 15,000 times 
units, and compute the Gini coefficient from the charging time of vehicles that have charged fully during 
the simulation. To reduce the effect of a transient regime, we consider only vehicles that are fully charged 
after iteration 1000. Vertical lines show the value of the critical points for max-flow (MF) and proportional 
fairness (PF) identified from the susceptibility xW for both networks. Symbols show average values over 
an ensemble of 25 runs and shaded areas represent 95% confidence intervals. 


that each vehicle gets a positive power allocation, thus the lower limit voltage constraint is satisfied 
at equality on the occupied node that is the most distant from the root on p. In max-flow, however, 
to maximise the aggregate power allocated to vehicles that can take all instantaneous power they 
are allocated (elastic demand), on a network with bounded voltage drops (i.e. capacity), implies 
also minimizing the power losses, and this is achieved by allocating all power on p to the closest 
occupied node from the root on that path. For max-flow, this implies vehicles on the path p further 
away from the root than the closest occupied node will only receive power after all vehicles on 
this node have left the network fully charged. In other words, under max-flow, users experience 
a charging time that depends strongly on their location on the network: vehicles close to the root 
charge faster, and vehicles on the tree leaves may take a very long time to charge. In contrast, 
under proportional fairness, the charging times are more homogeneous, because vehicles receive 
instantaneous powers that are also more uniform. 
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To characterise inequalities in the user experience, we analyse the Gini coelRcient of charging 


time. Origina 


defined as 11441] : 


ly devised as a measure of inequality in income distributions, the Gini coefficient is 


-I -I ^oo ^oo 

G = — E[|m- v|] = — I I \u-v\f(u)f(v)dudv, (11) 

2a/ 2ju Jo Jo 

where u and v are independent identically distributed random variables with probability density 
/ and mean /i. In other words, the Gini coefficient is one half of the mean difference in units 
of the mean. The difference between the two variables receives a small weight in the tail of the 
distribution, where f(u)f(v) is small, but a relatively large weight near the mode. Hence, G is 
more sensitive to changes near the mode than to changes in the tails. For a random sample (x,, 
i = 1,2,..., n), the empirical Gini coefficient, G, may be estimated by a sample mean 


G = 





( 12 ) 


The Gini coefficient is used as a measure of inequality, because a sample where the only non-zero 
value is x has /u = xin and hence G = (n - l)/n ^ 1 as n ^ oo, whereas G = 0 when all data 
points have the same value. 

We observe in Fig. |5] that the Gini coefficient of the charging time is larger in max-flow than in 
proportional fairness, for each of the networks. Moreover, the Gini coefficient increases faster in 
max-flow than in proportional fairness in the non-congested regime, showing that, when the system 
is stable, vehicles will experience a faster increase in the inequality of charging times in max-flow 
than in proportional fairness, with the increase of the vehicle arrival rate A. For comparison with 
well-known measures of income inequality, Sweden has a Gini of 0.26, the United States has 
a Gini of 0.41 and the Seychelles has the highest Gini of 0.66 1451]. The proportional fairness 
algorithm reaches a maximum Gini of 0.45, which is comparable with the level of inequality in 
the US society, and thus may be judged sociable acceptable. The max-flow algorithm, however, 
reaches a Gini of 0.91, which measures a level of inequality considerably higher than present in 
any contemporary society. 


IV. DISCUSSION 

In conclusion, we modelled the max-flow and proportional fairness protocols for the control 
of congestion caused by a fleet of vehicles charging on distribution networks. We analysed the 
second order phase transition that occurs with the increase of the number of electric vehicles that 
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arrive at the network with empty batteries to be charged, and found that the critical arrival rate 
Ac depends on the congestion control method. Indeed, we showed numerically on the 47-bus bus 
network that the onset of congestion takes place for larger values of A in proportional fairness than 
in max-flow. This result is surprising, because one would expect that, for a chosen arrival rate 
A, the maximisation of the aggregate instantaneous power would also lead to a maximisation of 
the energy (the time integral of power), and hence to a maximisation of the number of charged 
vehicles. This discovery, illustrates how the greediness of max-flow can be sub-optimal in relation 
to proportional fairness, which is an example of a fair allocation of instantaneous power. 

We analysed the inequality in the charging times as the vehicle arrival rate increases, and 
showed that charging times are considerably more equitable in proportional fairness than in max- 
flow. Indeed, vehicles close to the root get all the power allocation in max-flow, leaving other 
vehicles excluded from the network and unable to charge. Hence, proportional fairness is prefer¬ 
able to max-flow, not only because it does not exclude users from the network, but also because 
the charging times are more equitable, and it can serve a higher number of vehicles. In conclusion, 
proportional fairness is a promising candidate protocol to manage congestion in the charging of 
electric vehicles. 
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Appendix A: Voltage drop on one edge 


The angle 9 between V, and Vj is small in distribution networks 11361] (see Fig. 0, and hence 
the phases of V, and Vj are approximately the same, and can be chosen so the phasors have zero 
imaginary components. Since the phasors are real, we can derive the voltage drop from Kirchhoff’s 
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voltage law applied to the circuit in Fig.[IJb), 
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where the superscript asterisk denotes the complex conjugate transpose. 





FIG. 6. The difference lijZij between the y and Vj phasors, decomposed along the Vj vector and its 
orthogonal direction. The phase angle 6 difference between y and Vj is small, and hence the voltage drop 
can be approximated by Ay^ ^ % {jij Zi^. 


Appendix B: Active and reactive loads on a subtree 

From Kirchhoff’s current law, the active and reactive power consumed by the loads in the 
subtree rooted in node k can be computed as: 

N(f) 

Pm = ^iliOPliO + E E Pij(t), (Bl) 


and 


Q(\]ik) — QijiO^ 


(B2) 
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where is the active and the reactive power dissipated on a cable connecting nodes i 
and j. The complex power is given by: 


Sijit) = Pij(t) + iQijit) = iViit) - Vjit))!* = 
(Vi(t) - Vjity 

= {v,(t) - Vj(,)){v,(i) - v/t))' 


= {WM - W,j(t) - Wj:0) + W„(t)) 
Since, the voltages are real, Wij{t) = Wji(t), and thus 

Pij(t) = [Wn{t) - 2Wij{t) + Wjjit)) 
and 

Qij(t) = (Wuit)-2Wijit) + Wjj(t)) 




Rjj 

4 + 4’ 


X,- 


Rjj + xf: 


(B3) 
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(B5) 


Appendix C: Aggregation of vehicles at the nodes 


In proportional fairness, we maximise the sum of the logarithm of the instantaneous power 


allocated to electric vehicles: 

N{t) 


N(t) 


Ziog^'(^) = ZE A,7(0 log ■ 


Pi(0 


zrif A,(o 


l=l 


ie^+ 1=1 


(Cl) 


where Pi(t) is the instantaneous power allocated to electric vehicle /, and Pi the instantaneous 
power allocated to node i. To maximise Eq. (ICll) . we solve a problem with gradient and Hessian 


matrices that grow in size with the number of electric vehicles on the network. A more elRcient 
way to approach the problem is to aggregate cars for each node i, then solve the optimization 


problem for the nodes (as if they were ‘super-cars’), and finally distribute the power allocated to 


each node among the cars on the node. To do this, we remove constant terms in the objective 
function Eq. (ICll) . yielding: 

N{t) 

U{t) = ZE A,7(0 log Pi(0 = Z (C2) 
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